rna_data=read.csv("count.csv",header = TRUE,sep=" ",row.names = 1)
head(rna_data[1:5,1:5])
log_rna_data <- log2(rna_data + 1)

pdf("RNAseq_Histograms.pdf", width = 8, height = 6)
par(mfrow = c(2, 2))
for (i in 1:ncol(log_rna_data)) {
hist(log_rna_data[,i],                
       breaks = 50,                     # 分箱数量
       main = colnames(log_rna_data)[i], # 标题为样本名
       xlab = "Log2(Expression + 1)",
       ylab = "Frequency",
       col = "skyblue",
       border = "white",
       xlim = c(0, 20)
)
}
      
dev.off()  # 关闭图形设备

# 1. 提取单列数据
sample_data <- rna_data[, "PTB008D"]

# 2. 对数转换
log_sample_data <- log2(sample_data + 1)

# 3. 绘制直方图
hist(log_sample_data,
     breaks = 100,
     main = "PTB008D",
     xlab = "Log2($PTB008D + 1)",
     ylab = "Frequency",
     col = "lightblue",
     border = "navy"
)
dev.new()

data=read.csv("count.csv",header = TRUE,sep=" ",row.names = 1)#读取文件
pca <- prcomp(t(data))#执行PCA，需转置数据因为默认对行（样本）分析
x <- data.frame(t(data))#转换矩阵
plot(x$O00451,x$O00161)
library(ggplot2)
library(ggrepel)#防止标签重叠
x1 <- data.frame(pca$x)#包含主成分得分（样本在新坐标系的坐标）
x1$label <- rownames(x1)
ggplot(x1,aes(PC1,PC2)) +
  geom_point() +
  geom_label_repel(aes(label = label))
plot(pca$x[,1],pca$x[,2])
plot(pca)
biplot(pca,var.axes=FALSE)#biplot() 同时显示样本和变量（基因）关系

d = dist(t(data))#计算距离矩阵（默认欧氏距离）
clu <- hclust(d)#执行层次聚类（默认complete方法）
plot(clu)



load("volcano.RData")#加载包含差异分析结果的RData文件（通常包含基因表达和p值）
plotdata <- prostat# 将数据赋值给plotdata
colnames(plotdata)[2] <- c("log2fc")# 将第二列列名改为"log2fc"
plotdata$group <- "nosig"# 初始化分组为"nosig"（无显著差异）
pos = plotdata$log2fc>0.58 & plotdata$P<0.05# 标记上调基因（log2fc > 0.58且p < 0.05）
plotdata$group[pos] <- "up"
pos = plotdata$log2fc< -0.58 & plotdata$P<0.05# 标记下调基因（log2fc < -0.58且p < 0.05）
plotdata$group[pos] <- "down"
plotdata$label <- plotdata$ID# 添加标签列，仅保留差异基因的ID（非显著基因标签为空）
plotdata$label[plotdata$group == "nosig"] <- ""
library(ggplot2)
ggplot(plotdata,aes(x=log2fc,y=-log10(P))) +
  geom_point(aes(color = group)) +# 按分组着色
  scale_color_manual(values = c("green","grey","red"),# 颜色映射：下调=绿，无差异=灰，上调=红
                     limits = c("down","nosig","up")) +# 指定颜色顺序
  geom_label_repel(aes(label = label))# 添加防重叠的基因标签（仅差异基因）

x <- plotdata;
x$log10p <- -log10(x$P)计算-log10(p值)方便绘图
x$color <- x$group# 将分组转换为颜色代码
x$color[x$color == "nosig"] <- "grey"
x$color[x$color == "up"] <- "red"
x$color[x$color == "down"] <- "green"
plot(x$log2fc, x$log10p, type = "n",# 初始化绘图区域（type = "n"表示不绘制数据点）
     main="volcano plot",xlab="FC",ylab="-log10 p-value")
points(x$log2fc,x$log10p,col=x$color,pch=16,cex=1)
#points(1.5,6,col="black")
abline(v=log2(1.5),lty=3)# 添加阈值线
abline(v=log2(1/1.5),lty=3)
abline(h=-log10(0.05),lty=3)
legend(-1.5,1,c("up","nosig","down"),col = c("red","grey","green"),text.col ="black",# 添加图例
       pch = c(16,16,16),cex = c(0.4,0.4,0.4))
up <- subset(x,color == "red")# 为上调基因添加标签（y轴位置+0.2防止重叠）
text(up$log2fc,up$log10p+0.2,labels = up$ID)



